# Buget 0.1# Read the RDS file from the Data folderdata <-readRDS(here::here("Data/EDGE_full_results_averaged_budget0.1_replicates100.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Transform the data to sf object and projectdata_sf <-st_as_sf(data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the outline from the projected world mapoutline <-st_union(world_projected) %>%st_boundary()# Create the globe bounding boxglobe_bbox <-rbind(c(-180, -90), c(-180, 90), c(180, 90), c(180, -90), c(-180, -90))# Create the globe borderglobe_border <-st_polygon(list(globe_bbox)) %>%st_sfc(crs =4326) %>%st_sf(data.frame(rgn ='globe', geom = .)) %>% smoothr::densify(max_distance =0.5) %>%st_transform(crs = mcbryde_thomas_2)# Create the plotggplot() +# Add points for each planning unit, colored by prioritygeom_sf(data = data_sf, aes(color = Priority), size =0.5, alpha =0.7) +# Add the world map with a light gray fillgeom_sf(data = world_projected, fill ="lightgrey", color ="lightgrey", size =0.1) +# Add the black outline around the globegeom_sf(data = globe_border, fill =NA, color ="black", size =0.5) +# Add the globe border# Use a white to yellow to blue color gradientscale_color_gradientn(colors =c("white", "yellow", "darkblue"),values =c(0, 0.5, 1),name ="Priority",guide =guide_colorbar(barwidth =20, barheight =0.5, title.position ="top", title.hjust =0.5) ) +# Add labels and titlelabs(title ="Global Distribution of Conservation Priorities",subtitle ="Index: EDGE2, Budget: 0.1, Replicates: 100",x ="Longitude", y ="Latitude") +# Adjust themetheme_minimal() +theme(legend.position ="bottom",legend.direction ="horizontal",legend.box ="vertical",legend.margin =margin(t =20, r =0, b =0, l =0),legend.title =element_text(margin =margin(b =10)),panel.background =element_rect(fill ="white", color =NA),plot.background =element_rect(fill ="white", color =NA),panel.grid =element_blank() )
Code
# Save the plot#ggsave("priority_world_map_mcbryde_outline_EDGE2_0.1_100reps.png", width = 15, height = 10, dpi = 300)# Protection fraction summary# Load required librarieslibrary(ggplot2)library(jsonlite)library(here)library(gridExtra)# Read the dataprot_frac <-readRDS(here::here("Data/protect_fraction_summary_EDGE_0.1.rds"))sp <-fromJSON(here("Data", "shark_conservation_metrics_no_freshwater.json"))# Extract Species and EDGE2 from spSpecies <- sp$EDGE2$info$SpeciesEDGE2 <- sp$EDGE2$info$EDGE2# Add Species and EDGE2 to prot_fracprot_frac$Species <- Speciesprot_frac$EDGE2 <-as.numeric(EDGE2)# Create histogram for Mean_Protect_Fractionhist_protect <-ggplot(prot_frac, aes(x = Mean_Protect_Fraction)) +geom_histogram(binwidth =0.05, fill ="skyblue", color ="black") +theme_minimal() +labs(title ="Histogram of Mean Protect Fraction",x ="Mean Protect Fraction",y ="Count")# Create histogram for EDGE2hist_edge2 <-ggplot(prot_frac, aes(x = EDGE2)) +geom_histogram(binwidth =0.05, fill ="lightgreen", color ="black") +theme_minimal() +labs(title ="Histogram of EDGE2 Scores",x ="EDGE2 Score",y ="Count")# Create scatterplotscatter_plot <-ggplot(prot_frac, aes(x = EDGE2, y = Mean_Protect_Fraction)) +geom_point(alpha =0.6, color ="darkblue") +theme_minimal() +labs(title ="Scatterplot: EDGE2 vs Mean Protect Fraction",x ="EDGE2 Score",y ="Mean Protect Fraction")# Arrange plots in a gridgrid_plot <-grid.arrange( hist_protect, hist_edge2, scatter_plot,layout_matrix =rbind(c(1,2), c(3,3)),widths =c(1, 1),heights =c(1, 1))
Code
# Optionally, you can save the grid plot# ggsave("grid_plot.png", grid_plot, width = 12, height = 10)
FUSE
Code
# Read the RDS file from the Data folderdata <-readRDS(here::here("Data/FUSE_full_results_averaged_budget0.1_replicates100.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Transform the data to sf object and projectdata_sf <-st_as_sf(data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the outline from the projected world mapoutline <-st_union(world_projected) %>%st_boundary()# Create the plotggplot() +# Add points for each planning unit, colored by prioritygeom_sf(data = data_sf, aes(color = Priority), size =0.5, alpha =0.7) +# Add the world map with a light gray fillgeom_sf(data = world_projected, fill ="lightgrey", color ="lightgrey", size =0.1) +# Add the black outline around the globegeom_sf(data = globe_border, fill =NA, color ="black", size =0.5) +# Add the globe border# Use a white to yellow to blue color gradientscale_color_gradientn(colors =c("white", "yellow", "darkblue"),values =c(0, 0.5, 1),name ="Priority",guide =guide_colorbar(barwidth =20, barheight =0.5, title.position ="top", title.hjust =0.5) ) +# Add labels and titlelabs(title ="Global Distribution of Conservation Priorities",subtitle ="Index: FUSE, Budget: 0.1, Replicates: 100",x =NULL, y =NULL) +# Adjust themetheme_minimal() +theme(legend.position ="bottom",legend.direction ="horizontal",legend.box ="vertical",legend.margin =margin(t =20, r =0, b =0, l =0),legend.title =element_text(margin =margin(b =10)),panel.background =element_rect(fill ="white", color =NA),plot.background =element_rect(fill ="white", color =NA),panel.grid =element_blank() )
Code
# Protection fraction summary# Load required librarieslibrary(ggplot2)library(jsonlite)library(here)library(gridExtra)# Read the dataprot_frac <-readRDS(here::here("Data/protect_fraction_summary_FUSE_01.rds"))sp <-fromJSON(here("Data", "shark_conservation_metrics_no_freshwater.json"))# Extract Species and EDGE2 from spSpecies <- sp$FUSE$info$SpeciesFUSE <- sp$FUSE$info$FUSE# Add Species and EDGE2 to prot_fracprot_frac$Species <- Speciesprot_frac$FUSE <-as.numeric(FUSE)# Create histogram for Mean_Protect_Fractionhist_protect <-ggplot(prot_frac, aes(x = Mean_Protect_Fraction)) +geom_histogram(binwidth =0.05, fill ="skyblue", color ="black") +theme_minimal() +labs(title ="Histogram of Mean Protect Fraction",x ="Mean Protect Fraction",y ="Count")# Create histogram for EDGE2hist_fuse <-ggplot(prot_frac, aes(x = FUSE)) +geom_histogram(binwidth =0.05, fill ="lightgreen", color ="black") +theme_minimal() +labs(title ="Histogram of FUSE Scores",x ="FUSE Score",y ="Count")# Create scatterplotscatter_plot <-ggplot(prot_frac, aes(x = FUSE, y = Mean_Protect_Fraction)) +geom_point(alpha =0.6, color ="darkblue") +theme_minimal() +labs(title ="Scatterplot: FUSE vs Mean Protect Fraction",x ="FUSE Score",y ="Mean Protect Fraction")# Arrange plots in a gridgrid_plot <-grid.arrange( hist_protect, hist_fuse, scatter_plot,layout_matrix =rbind(c(1,2), c(3,3)),widths =c(1, 1),heights =c(1, 1))
EDGE2 and FUSE
Code
# RGB map ----# Read the RDS files from the Data folderedge_data <-readRDS(here::here("Data/EDGE_full_results_averaged_budget0.1_replicates100.rds"))fuse_data <-readRDS(here::here("Data/FUSE_full_results_averaged_budget0.1_replicates100.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Combine EDGE and FUSE datacombined_data <- edge_data %>%rename(EDGE_Priority = Priority) %>%left_join(fuse_data %>%rename(FUSE_Priority = Priority),by =c("Longitude", "Latitude"))# Normalize priorities to 0-1 rangecombined_data <- combined_data %>%mutate(EDGE_Priority_Norm = (EDGE_Priority -min(EDGE_Priority)) / (max(EDGE_Priority) -min(EDGE_Priority)),FUSE_Priority_Norm = (FUSE_Priority -min(FUSE_Priority)) / (max(FUSE_Priority) -min(FUSE_Priority)) )# Transform the data to sf object and projectdata_sf <-st_as_sf(combined_data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the color palettemap_pal_raw <-bi_pal(pal ='PurpleOr', dim =4, preview =FALSE)map_pal_mtx <-matrix(map_pal_raw, nrow =4, ncol =4)map_pal_mtx[3, ] <- colorspace::lighten(map_pal_mtx[3, ], .1)map_pal_mtx[2, ] <- colorspace::lighten(map_pal_mtx[2, ], .2)map_pal_mtx[1, ] <- colorspace::lighten(map_pal_mtx[1, ], .3)map_pal_mtx[ , 3] <- colorspace::lighten(map_pal_mtx[ , 3], .1)map_pal_mtx[ , 2] <- colorspace::lighten(map_pal_mtx[ , 2], .2)map_pal_mtx[ , 1] <- colorspace::lighten(map_pal_mtx[ , 1], .3)map_pal_mtx[1, 1] <-'#ffffee'map_pal <-as.vector(map_pal_mtx) %>%setNames(names(map_pal_raw))# Create a custom color mapping functionget_color <-function(edge, fuse) { edge_class <-cut(edge, breaks =c(-Inf, 0.25, 0.5, 0.75, Inf), labels =1:4) fuse_class <-cut(fuse, breaks =c(-Inf, 0.25, 0.5, 0.75, Inf), labels =1:4)return(map_pal[(as.numeric(fuse_class)-1)*4+as.numeric(edge_class)])}# Apply the function to get colorsdata_sf$new_color <-mapply(get_color, data_sf$EDGE_Priority_Norm, data_sf$FUSE_Priority_Norm)# Create the main plotmain_plot <-ggplot() +geom_sf(data = data_sf, aes(color = new_color), size =0.1, alpha =1) +geom_sf(data = world_projected, fill ="lightgray", color ="lightgray") +geom_sf(data = globe_border, fill =NA, color ="grey70", size =0.5) +# Add the globe borderscale_color_identity() +coord_sf() +theme_minimal() +labs(title ="Bivariate Map of EDGE and FUSE Priorities",x ="Longitude", y ="Latitude") +theme(plot.title =element_text(hjust =0.5))# Create the legendlegend_plot <-bi_legend(pal = map_pal, dim =4,xlab ='EDGE2',ylab ='FUSE')# Define the layout matrixlayout <-rbind(c(1, 1),c(2, 3))# Combine the main plot and legend using grid.arrangecombined_plot <-grid.arrange( main_plot, legend_plot,rectGrob(gp =gpar(col ="white")), # Blank spacelayout_matrix = layout,widths =c(0.25, 0.75), # Legend takes 15% width, main plot and blank space take 85%heights =c(0.8, 0.2) # Main plot takes 80% of height, legend row takes 20%)
Code
# Display the combined plot#print(combined_plot)# Save the combined plot#ggsave("bivariate_priority_map_with_legend_100reps_budget01.png", combined_plot, width = 15, height = 12, dpi = 300)
Budget 0.3
EDGE2
Code
# Buget 0.3# Read the RDS file from the Data folderdata <-readRDS(here::here("Data/EDGE_full_results_averaged_budget0.3_replicates30.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Transform the data to sf object and projectdata_sf <-st_as_sf(data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the outline from the projected world mapoutline <-st_union(world_projected) %>%st_boundary()# Create the globe bounding boxglobe_bbox <-rbind(c(-180, -90), c(-180, 90), c(180, 90), c(180, -90), c(-180, -90))# Create the globe borderglobe_border <-st_polygon(list(globe_bbox)) %>%st_sfc(crs =4326) %>%st_sf(data.frame(rgn ='globe', geom = .)) %>% smoothr::densify(max_distance =0.5) %>%st_transform(crs = mcbryde_thomas_2)# Create the plotggplot() +# Add points for each planning unit, colored by prioritygeom_sf(data = data_sf, aes(color = Priority), size =0.5, alpha =0.7) +# Add the world map with a light gray fillgeom_sf(data = world_projected, fill ="lightgrey", color ="lightgrey", size =0.1) +# Add the black outline around the globegeom_sf(data = globe_border, fill =NA, color ="black", size =0.5) +# Add the globe border# Use a white to yellow to blue color gradientscale_color_gradientn(colors =c("white", "yellow", "darkblue"),values =c(0, 0.5, 1),name ="Priority",guide =guide_colorbar(barwidth =20, barheight =0.5, title.position ="top", title.hjust =0.5) ) +# Add labels and titlelabs(title ="Global Distribution of Conservation Priorities",subtitle ="Index: EDGE2, Budget: 0.3 Replicates: 30",x ="Longitude", y ="Latitude") +# Adjust themetheme_minimal() +theme(legend.position ="bottom",legend.direction ="horizontal",legend.box ="vertical",legend.margin =margin(t =20, r =0, b =0, l =0),legend.title =element_text(margin =margin(b =10)),panel.background =element_rect(fill ="white", color =NA),plot.background =element_rect(fill ="white", color =NA),panel.grid =element_blank() )
Code
# Save the plot#ggsave("priority_world_map_mcbryde_outline_EDGE2_0.3_30reps.png", width = 15, height = 10, dpi = 300)# Protection fraction summary# Load required librarieslibrary(ggplot2)library(jsonlite)library(here)library(gridExtra)# Read the dataprot_frac <-readRDS(here::here("Data/protect_fraction_summary_EDGE_03.rds"))sp <-fromJSON(here("Data", "shark_conservation_metrics_no_freshwater.json"))# Extract Species and EDGE2 from spSpecies <- sp$EDGE2$info$SpeciesEDGE2 <- sp$EDGE2$info$EDGE2# Add Species and EDGE2 to prot_fracprot_frac$Species <- Speciesprot_frac$EDGE2 <-as.numeric(EDGE2)# Create histogram for Mean_Protect_Fractionhist_protect <-ggplot(prot_frac, aes(x = Mean_Protect_Fraction)) +geom_histogram(binwidth =0.05, fill ="skyblue", color ="black") +theme_minimal() +labs(title ="Histogram of Mean Protect Fraction",x ="Mean Protect Fraction",y ="Count")# Create histogram for EDGE2hist_edge2 <-ggplot(prot_frac, aes(x = EDGE2)) +geom_histogram(binwidth =0.05, fill ="lightgreen", color ="black") +theme_minimal() +labs(title ="Histogram of EDGE2 Scores",x ="EDGE2 Score",y ="Count")# Create scatterplotscatter_plot <-ggplot(prot_frac, aes(x = EDGE2, y = Mean_Protect_Fraction)) +geom_point(alpha =0.6, color ="darkblue") +theme_minimal() +labs(title ="Scatterplot: EDGE2 vs Mean Protect Fraction",x ="EDGE2 Score",y ="Mean Protect Fraction")# Arrange plots in a gridgrid_plot <-grid.arrange( hist_protect, hist_edge2, scatter_plot,layout_matrix =rbind(c(1,2), c(3,3)),widths =c(1, 1),heights =c(1, 1))
FUSE
Code
# Buget 0.3# Read the RDS file from the Data folderdata <-readRDS(here::here("Data/FUSE_full_results_averaged_budget0.3_replicates56.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Transform the data to sf object and projectdata_sf <-st_as_sf(data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the outline from the projected world mapoutline <-st_union(world_projected) %>%st_boundary()# Create the globe bounding boxglobe_bbox <-rbind(c(-180, -90), c(-180, 90), c(180, 90), c(180, -90), c(-180, -90))# Create the globe borderglobe_border <-st_polygon(list(globe_bbox)) %>%st_sfc(crs =4326) %>%st_sf(data.frame(rgn ='globe', geom = .)) %>% smoothr::densify(max_distance =0.5) %>%st_transform(crs = mcbryde_thomas_2)# Create the plotggplot() +# Add points for each planning unit, colored by prioritygeom_sf(data = data_sf, aes(color = Priority), size =0.5, alpha =0.7) +# Add the world map with a light gray fillgeom_sf(data = world_projected, fill ="lightgrey", color ="lightgrey", size =0.1) +# Add the black outline around the globegeom_sf(data = globe_border, fill =NA, color ="black", size =0.5) +# Add the globe border# Use a white to yellow to blue color gradientscale_color_gradientn(colors =c("white", "yellow", "darkblue"),values =c(0, 0.5, 1),name ="Priority",guide =guide_colorbar(barwidth =20, barheight =0.5, title.position ="top", title.hjust =0.5) ) +# Add labels and titlelabs(title ="Global Distribution of Conservation Priorities",subtitle ="Index: FUSE, Budget: 0.3 Replicates: 30",x ="Longitude", y ="Latitude") +# Adjust themetheme_minimal() +theme(legend.position ="bottom",legend.direction ="horizontal",legend.box ="vertical",legend.margin =margin(t =20, r =0, b =0, l =0),legend.title =element_text(margin =margin(b =10)),panel.background =element_rect(fill ="white", color =NA),plot.background =element_rect(fill ="white", color =NA),panel.grid =element_blank() )
Code
# Save the plot#ggsave("priority_world_map_mcbryde_outline_FUSE_0.3_30reps.png", width = 15, height = 10, dpi = 300)# Protection fraction summary# Load required librarieslibrary(ggplot2)library(jsonlite)library(here)library(gridExtra)# Read the dataprot_frac <-readRDS(here::here("Data/protect_fraction_summary_EDGE_03.rds"))sp <-fromJSON(here("Data", "shark_conservation_metrics_no_freshwater.json"))# Extract Species and FUSE from spSpecies <- sp$FUSE$info$SpeciesFUSE <- sp$FUSE$info$FUSE# Add Species and FUSE to prot_fracprot_frac$Species <- Speciesprot_frac$FUSE <-as.numeric(FUSE)# Create histogram for Mean_Protect_Fractionhist_protect <-ggplot(prot_frac, aes(x = Mean_Protect_Fraction)) +geom_histogram(binwidth =0.05, fill ="skyblue", color ="black") +theme_minimal() +labs(title ="Histogram of Mean Protect Fraction",x ="Mean Protect Fraction",y ="Count")# Create histogram for FUSEhist_FUSE <-ggplot(prot_frac, aes(x = FUSE)) +geom_histogram(binwidth =0.05, fill ="lightgreen", color ="black") +theme_minimal() +labs(title ="Histogram of FUSE Scores",x ="FUSE Score",y ="Count")# Create scatterplotscatter_plot <-ggplot(prot_frac, aes(x = FUSE, y = Mean_Protect_Fraction)) +geom_point(alpha =0.6, color ="darkblue") +theme_minimal() +labs(title ="Scatterplot: FUSE vs Mean Protect Fraction",x ="FUSE Score",y ="Mean Protect Fraction")# Arrange plots in a gridgrid_plot <-grid.arrange( hist_protect, hist_FUSE, scatter_plot,layout_matrix =rbind(c(1,2), c(3,3)),widths =c(1, 1),heights =c(1, 1))
EDGE2 and FUSE
Code
# RGB map ----# Read the RDS files from the Data folderedge_data <-readRDS(here::here("Data/EDGE_full_results_averaged_budget0.3_replicates30.rds"))fuse_data <-readRDS(here::here("Data/FUSE_full_results_averaged_budget0.3_replicates56.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Combine EDGE and FUSE datacombined_data <- edge_data %>%rename(EDGE_Priority = Priority) %>%left_join(fuse_data %>%rename(FUSE_Priority = Priority),by =c("Longitude", "Latitude"))# Normalize priorities to 0-1 rangecombined_data <- combined_data %>%mutate(EDGE_Priority_Norm = (EDGE_Priority -min(EDGE_Priority)) / (max(EDGE_Priority) -min(EDGE_Priority)),FUSE_Priority_Norm = (FUSE_Priority -min(FUSE_Priority)) / (max(FUSE_Priority) -min(FUSE_Priority)) )# Transform the data to sf object and projectdata_sf <-st_as_sf(combined_data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the color palettemap_pal_raw <-bi_pal(pal ='PurpleOr', dim =4, preview =FALSE)map_pal_mtx <-matrix(map_pal_raw, nrow =4, ncol =4)map_pal_mtx[3, ] <- colorspace::lighten(map_pal_mtx[3, ], .1)map_pal_mtx[2, ] <- colorspace::lighten(map_pal_mtx[2, ], .2)map_pal_mtx[1, ] <- colorspace::lighten(map_pal_mtx[1, ], .3)map_pal_mtx[ , 3] <- colorspace::lighten(map_pal_mtx[ , 3], .1)map_pal_mtx[ , 2] <- colorspace::lighten(map_pal_mtx[ , 2], .2)map_pal_mtx[ , 1] <- colorspace::lighten(map_pal_mtx[ , 1], .3)map_pal_mtx[1, 1] <-'#ffffee'map_pal <-as.vector(map_pal_mtx) %>%setNames(names(map_pal_raw))# Create a custom color mapping functionget_color <-function(edge, fuse) { edge_class <-cut(edge, breaks =c(-Inf, 0.25, 0.5, 0.75, Inf), labels =1:4) fuse_class <-cut(fuse, breaks =c(-Inf, 0.25, 0.5, 0.75, Inf), labels =1:4)return(map_pal[(as.numeric(fuse_class)-1)*4+as.numeric(edge_class)])}# Apply the function to get colorsdata_sf$new_color <-mapply(get_color, data_sf$EDGE_Priority_Norm, data_sf$FUSE_Priority_Norm)# Create the main plotmain_plot <-ggplot() +geom_sf(data = data_sf, aes(color = new_color), size =0.1, alpha =1) +geom_sf(data = world_projected, fill ="lightgray", color ="lightgray") +geom_sf(data = globe_border, fill =NA, color ="grey70", size =0.5) +# Add the globe borderscale_color_identity() +coord_sf() +theme_minimal() +labs(title ="Bivariate Map of EDGE and FUSE Priorities",x ="Longitude", y ="Latitude") +theme(plot.title =element_text(hjust =0.5))# Create the legendlegend_plot <-bi_legend(pal = map_pal, dim =4,xlab ='EDGE2',ylab ='FUSE')# Define the layout matrixlayout <-rbind(c(1, 1),c(2, 3))# Combine the main plot and legend using grid.arrangecombined_plot <-grid.arrange( main_plot, legend_plot,rectGrob(gp =gpar(col ="white")), # Blank spacelayout_matrix = layout,widths =c(0.25, 0.75), # Legend takes 15% width, main plot and blank space take 85%heights =c(0.8, 0.2) # Main plot takes 80% of height, legend row takes 20%)
Code
# Display the combined plot#print(combined_plot)# Save the combined plot#ggsave("bivariate_priority_map_with_legend_100reps_budget01.png", combined_plot, width = 15, height = 12, dpi = 300)
Source Code
---title: "CAPTAIN visualisation workflow"author: "Théophile L. Mouton"date: "September 24, 2024"format: html: toc: true toc-location: right css: custom.css output-file: "CAPTAIN_workflow.html" self-contained: true code-fold: true code-tools: trueeditor: visualexecute: warning: false message: false echo: true---# Visualise CAPTAIN results### R libraries```{r}library(readr)library(ggplot2)library(sf)library(rnaturalearth)library(rnaturalearthdata)library(dplyr)library(gridExtra)library(biscale)library(colorspace)library(grid)```## Budget 0.1 ### EDGE2```{r}# Buget 0.1# Read the RDS file from the Data folderdata <-readRDS(here::here("Data/EDGE_full_results_averaged_budget0.1_replicates100.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Transform the data to sf object and projectdata_sf <-st_as_sf(data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the outline from the projected world mapoutline <-st_union(world_projected) %>%st_boundary()# Create the globe bounding boxglobe_bbox <-rbind(c(-180, -90), c(-180, 90), c(180, 90), c(180, -90), c(-180, -90))# Create the globe borderglobe_border <-st_polygon(list(globe_bbox)) %>%st_sfc(crs =4326) %>%st_sf(data.frame(rgn ='globe', geom = .)) %>% smoothr::densify(max_distance =0.5) %>%st_transform(crs = mcbryde_thomas_2)# Create the plotggplot() +# Add points for each planning unit, colored by prioritygeom_sf(data = data_sf, aes(color = Priority), size =0.5, alpha =0.7) +# Add the world map with a light gray fillgeom_sf(data = world_projected, fill ="lightgrey", color ="lightgrey", size =0.1) +# Add the black outline around the globegeom_sf(data = globe_border, fill =NA, color ="black", size =0.5) +# Add the globe border# Use a white to yellow to blue color gradientscale_color_gradientn(colors =c("white", "yellow", "darkblue"),values =c(0, 0.5, 1),name ="Priority",guide =guide_colorbar(barwidth =20, barheight =0.5, title.position ="top", title.hjust =0.5) ) +# Add labels and titlelabs(title ="Global Distribution of Conservation Priorities",subtitle ="Index: EDGE2, Budget: 0.1, Replicates: 100",x ="Longitude", y ="Latitude") +# Adjust themetheme_minimal() +theme(legend.position ="bottom",legend.direction ="horizontal",legend.box ="vertical",legend.margin =margin(t =20, r =0, b =0, l =0),legend.title =element_text(margin =margin(b =10)),panel.background =element_rect(fill ="white", color =NA),plot.background =element_rect(fill ="white", color =NA),panel.grid =element_blank() )# Save the plot#ggsave("priority_world_map_mcbryde_outline_EDGE2_0.1_100reps.png", width = 15, height = 10, dpi = 300)# Protection fraction summary# Load required librarieslibrary(ggplot2)library(jsonlite)library(here)library(gridExtra)# Read the dataprot_frac <-readRDS(here::here("Data/protect_fraction_summary_EDGE_0.1.rds"))sp <-fromJSON(here("Data", "shark_conservation_metrics_no_freshwater.json"))# Extract Species and EDGE2 from spSpecies <- sp$EDGE2$info$SpeciesEDGE2 <- sp$EDGE2$info$EDGE2# Add Species and EDGE2 to prot_fracprot_frac$Species <- Speciesprot_frac$EDGE2 <-as.numeric(EDGE2)# Create histogram for Mean_Protect_Fractionhist_protect <-ggplot(prot_frac, aes(x = Mean_Protect_Fraction)) +geom_histogram(binwidth =0.05, fill ="skyblue", color ="black") +theme_minimal() +labs(title ="Histogram of Mean Protect Fraction",x ="Mean Protect Fraction",y ="Count")# Create histogram for EDGE2hist_edge2 <-ggplot(prot_frac, aes(x = EDGE2)) +geom_histogram(binwidth =0.05, fill ="lightgreen", color ="black") +theme_minimal() +labs(title ="Histogram of EDGE2 Scores",x ="EDGE2 Score",y ="Count")# Create scatterplotscatter_plot <-ggplot(prot_frac, aes(x = EDGE2, y = Mean_Protect_Fraction)) +geom_point(alpha =0.6, color ="darkblue") +theme_minimal() +labs(title ="Scatterplot: EDGE2 vs Mean Protect Fraction",x ="EDGE2 Score",y ="Mean Protect Fraction")# Arrange plots in a gridgrid_plot <-grid.arrange( hist_protect, hist_edge2, scatter_plot,layout_matrix =rbind(c(1,2), c(3,3)),widths =c(1, 1),heights =c(1, 1))# Optionally, you can save the grid plot# ggsave("grid_plot.png", grid_plot, width = 12, height = 10)```### FUSE```{r}# Read the RDS file from the Data folderdata <-readRDS(here::here("Data/FUSE_full_results_averaged_budget0.1_replicates100.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Transform the data to sf object and projectdata_sf <-st_as_sf(data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the outline from the projected world mapoutline <-st_union(world_projected) %>%st_boundary()# Create the plotggplot() +# Add points for each planning unit, colored by prioritygeom_sf(data = data_sf, aes(color = Priority), size =0.5, alpha =0.7) +# Add the world map with a light gray fillgeom_sf(data = world_projected, fill ="lightgrey", color ="lightgrey", size =0.1) +# Add the black outline around the globegeom_sf(data = globe_border, fill =NA, color ="black", size =0.5) +# Add the globe border# Use a white to yellow to blue color gradientscale_color_gradientn(colors =c("white", "yellow", "darkblue"),values =c(0, 0.5, 1),name ="Priority",guide =guide_colorbar(barwidth =20, barheight =0.5, title.position ="top", title.hjust =0.5) ) +# Add labels and titlelabs(title ="Global Distribution of Conservation Priorities",subtitle ="Index: FUSE, Budget: 0.1, Replicates: 100",x =NULL, y =NULL) +# Adjust themetheme_minimal() +theme(legend.position ="bottom",legend.direction ="horizontal",legend.box ="vertical",legend.margin =margin(t =20, r =0, b =0, l =0),legend.title =element_text(margin =margin(b =10)),panel.background =element_rect(fill ="white", color =NA),plot.background =element_rect(fill ="white", color =NA),panel.grid =element_blank() )# Protection fraction summary# Load required librarieslibrary(ggplot2)library(jsonlite)library(here)library(gridExtra)# Read the dataprot_frac <-readRDS(here::here("Data/protect_fraction_summary_FUSE_01.rds"))sp <-fromJSON(here("Data", "shark_conservation_metrics_no_freshwater.json"))# Extract Species and EDGE2 from spSpecies <- sp$FUSE$info$SpeciesFUSE <- sp$FUSE$info$FUSE# Add Species and EDGE2 to prot_fracprot_frac$Species <- Speciesprot_frac$FUSE <-as.numeric(FUSE)# Create histogram for Mean_Protect_Fractionhist_protect <-ggplot(prot_frac, aes(x = Mean_Protect_Fraction)) +geom_histogram(binwidth =0.05, fill ="skyblue", color ="black") +theme_minimal() +labs(title ="Histogram of Mean Protect Fraction",x ="Mean Protect Fraction",y ="Count")# Create histogram for EDGE2hist_fuse <-ggplot(prot_frac, aes(x = FUSE)) +geom_histogram(binwidth =0.05, fill ="lightgreen", color ="black") +theme_minimal() +labs(title ="Histogram of FUSE Scores",x ="FUSE Score",y ="Count")# Create scatterplotscatter_plot <-ggplot(prot_frac, aes(x = FUSE, y = Mean_Protect_Fraction)) +geom_point(alpha =0.6, color ="darkblue") +theme_minimal() +labs(title ="Scatterplot: FUSE vs Mean Protect Fraction",x ="FUSE Score",y ="Mean Protect Fraction")# Arrange plots in a gridgrid_plot <-grid.arrange( hist_protect, hist_fuse, scatter_plot,layout_matrix =rbind(c(1,2), c(3,3)),widths =c(1, 1),heights =c(1, 1))```### EDGE2 and FUSE```{r}# RGB map ----# Read the RDS files from the Data folderedge_data <-readRDS(here::here("Data/EDGE_full_results_averaged_budget0.1_replicates100.rds"))fuse_data <-readRDS(here::here("Data/FUSE_full_results_averaged_budget0.1_replicates100.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Combine EDGE and FUSE datacombined_data <- edge_data %>%rename(EDGE_Priority = Priority) %>%left_join(fuse_data %>%rename(FUSE_Priority = Priority),by =c("Longitude", "Latitude"))# Normalize priorities to 0-1 rangecombined_data <- combined_data %>%mutate(EDGE_Priority_Norm = (EDGE_Priority -min(EDGE_Priority)) / (max(EDGE_Priority) -min(EDGE_Priority)),FUSE_Priority_Norm = (FUSE_Priority -min(FUSE_Priority)) / (max(FUSE_Priority) -min(FUSE_Priority)) )# Transform the data to sf object and projectdata_sf <-st_as_sf(combined_data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the color palettemap_pal_raw <-bi_pal(pal ='PurpleOr', dim =4, preview =FALSE)map_pal_mtx <-matrix(map_pal_raw, nrow =4, ncol =4)map_pal_mtx[3, ] <- colorspace::lighten(map_pal_mtx[3, ], .1)map_pal_mtx[2, ] <- colorspace::lighten(map_pal_mtx[2, ], .2)map_pal_mtx[1, ] <- colorspace::lighten(map_pal_mtx[1, ], .3)map_pal_mtx[ , 3] <- colorspace::lighten(map_pal_mtx[ , 3], .1)map_pal_mtx[ , 2] <- colorspace::lighten(map_pal_mtx[ , 2], .2)map_pal_mtx[ , 1] <- colorspace::lighten(map_pal_mtx[ , 1], .3)map_pal_mtx[1, 1] <-'#ffffee'map_pal <-as.vector(map_pal_mtx) %>%setNames(names(map_pal_raw))# Create a custom color mapping functionget_color <-function(edge, fuse) { edge_class <-cut(edge, breaks =c(-Inf, 0.25, 0.5, 0.75, Inf), labels =1:4) fuse_class <-cut(fuse, breaks =c(-Inf, 0.25, 0.5, 0.75, Inf), labels =1:4)return(map_pal[(as.numeric(fuse_class)-1)*4+as.numeric(edge_class)])}# Apply the function to get colorsdata_sf$new_color <-mapply(get_color, data_sf$EDGE_Priority_Norm, data_sf$FUSE_Priority_Norm)# Create the main plotmain_plot <-ggplot() +geom_sf(data = data_sf, aes(color = new_color), size =0.1, alpha =1) +geom_sf(data = world_projected, fill ="lightgray", color ="lightgray") +geom_sf(data = globe_border, fill =NA, color ="grey70", size =0.5) +# Add the globe borderscale_color_identity() +coord_sf() +theme_minimal() +labs(title ="Bivariate Map of EDGE and FUSE Priorities",x ="Longitude", y ="Latitude") +theme(plot.title =element_text(hjust =0.5))# Create the legendlegend_plot <-bi_legend(pal = map_pal, dim =4,xlab ='EDGE2',ylab ='FUSE')# Define the layout matrixlayout <-rbind(c(1, 1),c(2, 3))# Combine the main plot and legend using grid.arrangecombined_plot <-grid.arrange( main_plot, legend_plot,rectGrob(gp =gpar(col ="white")), # Blank spacelayout_matrix = layout,widths =c(0.25, 0.75), # Legend takes 15% width, main plot and blank space take 85%heights =c(0.8, 0.2) # Main plot takes 80% of height, legend row takes 20%)# Display the combined plot#print(combined_plot)# Save the combined plot#ggsave("bivariate_priority_map_with_legend_100reps_budget01.png", combined_plot, width = 15, height = 12, dpi = 300)```## Budget 0.3 ### EDGE2```{r}# Buget 0.3# Read the RDS file from the Data folderdata <-readRDS(here::here("Data/EDGE_full_results_averaged_budget0.3_replicates30.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Transform the data to sf object and projectdata_sf <-st_as_sf(data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the outline from the projected world mapoutline <-st_union(world_projected) %>%st_boundary()# Create the globe bounding boxglobe_bbox <-rbind(c(-180, -90), c(-180, 90), c(180, 90), c(180, -90), c(-180, -90))# Create the globe borderglobe_border <-st_polygon(list(globe_bbox)) %>%st_sfc(crs =4326) %>%st_sf(data.frame(rgn ='globe', geom = .)) %>% smoothr::densify(max_distance =0.5) %>%st_transform(crs = mcbryde_thomas_2)# Create the plotggplot() +# Add points for each planning unit, colored by prioritygeom_sf(data = data_sf, aes(color = Priority), size =0.5, alpha =0.7) +# Add the world map with a light gray fillgeom_sf(data = world_projected, fill ="lightgrey", color ="lightgrey", size =0.1) +# Add the black outline around the globegeom_sf(data = globe_border, fill =NA, color ="black", size =0.5) +# Add the globe border# Use a white to yellow to blue color gradientscale_color_gradientn(colors =c("white", "yellow", "darkblue"),values =c(0, 0.5, 1),name ="Priority",guide =guide_colorbar(barwidth =20, barheight =0.5, title.position ="top", title.hjust =0.5) ) +# Add labels and titlelabs(title ="Global Distribution of Conservation Priorities",subtitle ="Index: EDGE2, Budget: 0.3 Replicates: 30",x ="Longitude", y ="Latitude") +# Adjust themetheme_minimal() +theme(legend.position ="bottom",legend.direction ="horizontal",legend.box ="vertical",legend.margin =margin(t =20, r =0, b =0, l =0),legend.title =element_text(margin =margin(b =10)),panel.background =element_rect(fill ="white", color =NA),plot.background =element_rect(fill ="white", color =NA),panel.grid =element_blank() )# Save the plot#ggsave("priority_world_map_mcbryde_outline_EDGE2_0.3_30reps.png", width = 15, height = 10, dpi = 300)# Protection fraction summary# Load required librarieslibrary(ggplot2)library(jsonlite)library(here)library(gridExtra)# Read the dataprot_frac <-readRDS(here::here("Data/protect_fraction_summary_EDGE_03.rds"))sp <-fromJSON(here("Data", "shark_conservation_metrics_no_freshwater.json"))# Extract Species and EDGE2 from spSpecies <- sp$EDGE2$info$SpeciesEDGE2 <- sp$EDGE2$info$EDGE2# Add Species and EDGE2 to prot_fracprot_frac$Species <- Speciesprot_frac$EDGE2 <-as.numeric(EDGE2)# Create histogram for Mean_Protect_Fractionhist_protect <-ggplot(prot_frac, aes(x = Mean_Protect_Fraction)) +geom_histogram(binwidth =0.05, fill ="skyblue", color ="black") +theme_minimal() +labs(title ="Histogram of Mean Protect Fraction",x ="Mean Protect Fraction",y ="Count")# Create histogram for EDGE2hist_edge2 <-ggplot(prot_frac, aes(x = EDGE2)) +geom_histogram(binwidth =0.05, fill ="lightgreen", color ="black") +theme_minimal() +labs(title ="Histogram of EDGE2 Scores",x ="EDGE2 Score",y ="Count")# Create scatterplotscatter_plot <-ggplot(prot_frac, aes(x = EDGE2, y = Mean_Protect_Fraction)) +geom_point(alpha =0.6, color ="darkblue") +theme_minimal() +labs(title ="Scatterplot: EDGE2 vs Mean Protect Fraction",x ="EDGE2 Score",y ="Mean Protect Fraction")# Arrange plots in a gridgrid_plot <-grid.arrange( hist_protect, hist_edge2, scatter_plot,layout_matrix =rbind(c(1,2), c(3,3)),widths =c(1, 1),heights =c(1, 1))```### FUSE```{r}# Buget 0.3# Read the RDS file from the Data folderdata <-readRDS(here::here("Data/FUSE_full_results_averaged_budget0.3_replicates56.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Transform the data to sf object and projectdata_sf <-st_as_sf(data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the outline from the projected world mapoutline <-st_union(world_projected) %>%st_boundary()# Create the globe bounding boxglobe_bbox <-rbind(c(-180, -90), c(-180, 90), c(180, 90), c(180, -90), c(-180, -90))# Create the globe borderglobe_border <-st_polygon(list(globe_bbox)) %>%st_sfc(crs =4326) %>%st_sf(data.frame(rgn ='globe', geom = .)) %>% smoothr::densify(max_distance =0.5) %>%st_transform(crs = mcbryde_thomas_2)# Create the plotggplot() +# Add points for each planning unit, colored by prioritygeom_sf(data = data_sf, aes(color = Priority), size =0.5, alpha =0.7) +# Add the world map with a light gray fillgeom_sf(data = world_projected, fill ="lightgrey", color ="lightgrey", size =0.1) +# Add the black outline around the globegeom_sf(data = globe_border, fill =NA, color ="black", size =0.5) +# Add the globe border# Use a white to yellow to blue color gradientscale_color_gradientn(colors =c("white", "yellow", "darkblue"),values =c(0, 0.5, 1),name ="Priority",guide =guide_colorbar(barwidth =20, barheight =0.5, title.position ="top", title.hjust =0.5) ) +# Add labels and titlelabs(title ="Global Distribution of Conservation Priorities",subtitle ="Index: FUSE, Budget: 0.3 Replicates: 30",x ="Longitude", y ="Latitude") +# Adjust themetheme_minimal() +theme(legend.position ="bottom",legend.direction ="horizontal",legend.box ="vertical",legend.margin =margin(t =20, r =0, b =0, l =0),legend.title =element_text(margin =margin(b =10)),panel.background =element_rect(fill ="white", color =NA),plot.background =element_rect(fill ="white", color =NA),panel.grid =element_blank() )# Save the plot#ggsave("priority_world_map_mcbryde_outline_FUSE_0.3_30reps.png", width = 15, height = 10, dpi = 300)# Protection fraction summary# Load required librarieslibrary(ggplot2)library(jsonlite)library(here)library(gridExtra)# Read the dataprot_frac <-readRDS(here::here("Data/protect_fraction_summary_EDGE_03.rds"))sp <-fromJSON(here("Data", "shark_conservation_metrics_no_freshwater.json"))# Extract Species and FUSE from spSpecies <- sp$FUSE$info$SpeciesFUSE <- sp$FUSE$info$FUSE# Add Species and FUSE to prot_fracprot_frac$Species <- Speciesprot_frac$FUSE <-as.numeric(FUSE)# Create histogram for Mean_Protect_Fractionhist_protect <-ggplot(prot_frac, aes(x = Mean_Protect_Fraction)) +geom_histogram(binwidth =0.05, fill ="skyblue", color ="black") +theme_minimal() +labs(title ="Histogram of Mean Protect Fraction",x ="Mean Protect Fraction",y ="Count")# Create histogram for FUSEhist_FUSE <-ggplot(prot_frac, aes(x = FUSE)) +geom_histogram(binwidth =0.05, fill ="lightgreen", color ="black") +theme_minimal() +labs(title ="Histogram of FUSE Scores",x ="FUSE Score",y ="Count")# Create scatterplotscatter_plot <-ggplot(prot_frac, aes(x = FUSE, y = Mean_Protect_Fraction)) +geom_point(alpha =0.6, color ="darkblue") +theme_minimal() +labs(title ="Scatterplot: FUSE vs Mean Protect Fraction",x ="FUSE Score",y ="Mean Protect Fraction")# Arrange plots in a gridgrid_plot <-grid.arrange( hist_protect, hist_FUSE, scatter_plot,layout_matrix =rbind(c(1,2), c(3,3)),widths =c(1, 1),heights =c(1, 1))```### EDGE2 and FUSE```{r}# RGB map ----# Read the RDS files from the Data folderedge_data <-readRDS(here::here("Data/EDGE_full_results_averaged_budget0.3_replicates30.rds"))fuse_data <-readRDS(here::here("Data/FUSE_full_results_averaged_budget0.3_replicates56.rds"))# Get world map dataworld <-ne_countries(scale ="medium", returnclass ="sf")# Define the McBryde-Thomas 2 projectionmcbryde_thomas_2 <-"+proj=mbt_s"# Combine EDGE and FUSE datacombined_data <- edge_data %>%rename(EDGE_Priority = Priority) %>%left_join(fuse_data %>%rename(FUSE_Priority = Priority),by =c("Longitude", "Latitude"))# Normalize priorities to 0-1 rangecombined_data <- combined_data %>%mutate(EDGE_Priority_Norm = (EDGE_Priority -min(EDGE_Priority)) / (max(EDGE_Priority) -min(EDGE_Priority)),FUSE_Priority_Norm = (FUSE_Priority -min(FUSE_Priority)) / (max(FUSE_Priority) -min(FUSE_Priority)) )# Transform the data to sf object and projectdata_sf <-st_as_sf(combined_data, coords =c("Longitude", "Latitude"), crs =4326) %>%st_transform(crs = mcbryde_thomas_2)# Project the world mapworld_projected <-st_transform(world, crs = mcbryde_thomas_2)# Create the color palettemap_pal_raw <-bi_pal(pal ='PurpleOr', dim =4, preview =FALSE)map_pal_mtx <-matrix(map_pal_raw, nrow =4, ncol =4)map_pal_mtx[3, ] <- colorspace::lighten(map_pal_mtx[3, ], .1)map_pal_mtx[2, ] <- colorspace::lighten(map_pal_mtx[2, ], .2)map_pal_mtx[1, ] <- colorspace::lighten(map_pal_mtx[1, ], .3)map_pal_mtx[ , 3] <- colorspace::lighten(map_pal_mtx[ , 3], .1)map_pal_mtx[ , 2] <- colorspace::lighten(map_pal_mtx[ , 2], .2)map_pal_mtx[ , 1] <- colorspace::lighten(map_pal_mtx[ , 1], .3)map_pal_mtx[1, 1] <-'#ffffee'map_pal <-as.vector(map_pal_mtx) %>%setNames(names(map_pal_raw))# Create a custom color mapping functionget_color <-function(edge, fuse) { edge_class <-cut(edge, breaks =c(-Inf, 0.25, 0.5, 0.75, Inf), labels =1:4) fuse_class <-cut(fuse, breaks =c(-Inf, 0.25, 0.5, 0.75, Inf), labels =1:4)return(map_pal[(as.numeric(fuse_class)-1)*4+as.numeric(edge_class)])}# Apply the function to get colorsdata_sf$new_color <-mapply(get_color, data_sf$EDGE_Priority_Norm, data_sf$FUSE_Priority_Norm)# Create the main plotmain_plot <-ggplot() +geom_sf(data = data_sf, aes(color = new_color), size =0.1, alpha =1) +geom_sf(data = world_projected, fill ="lightgray", color ="lightgray") +geom_sf(data = globe_border, fill =NA, color ="grey70", size =0.5) +# Add the globe borderscale_color_identity() +coord_sf() +theme_minimal() +labs(title ="Bivariate Map of EDGE and FUSE Priorities",x ="Longitude", y ="Latitude") +theme(plot.title =element_text(hjust =0.5))# Create the legendlegend_plot <-bi_legend(pal = map_pal, dim =4,xlab ='EDGE2',ylab ='FUSE')# Define the layout matrixlayout <-rbind(c(1, 1),c(2, 3))# Combine the main plot and legend using grid.arrangecombined_plot <-grid.arrange( main_plot, legend_plot,rectGrob(gp =gpar(col ="white")), # Blank spacelayout_matrix = layout,widths =c(0.25, 0.75), # Legend takes 15% width, main plot and blank space take 85%heights =c(0.8, 0.2) # Main plot takes 80% of height, legend row takes 20%)# Display the combined plot#print(combined_plot)# Save the combined plot#ggsave("bivariate_priority_map_with_legend_100reps_budget01.png", combined_plot, width = 15, height = 12, dpi = 300)```